//
// Created by 邹迪凯 on 2021/7/19.
//

#ifndef LWPT_UTIL_H
#define LWPT_UTIL_H

#include <cmath>
#include <iostream>
#include <vector>
#include <openssl/aes.h>
#include <openssl/evp.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>
#include <gsl/gsl_integration.h>
#define INF 1000000
//#include <boost/geometry.hpp>
//#pragma comment(lib, "./libidld.dylib" )
//#pragma comment(lib, "mclmcrrt.dylib")

//#include "libidld.h"

using namespace std;

const int PRIME_COUNT = 4;

#define PI 3.1415926535897932384626433832795
#define EARTH_RADIUS 6378.137 //地球半径 KM

template<class T>
void swap(vector<T> &vec, int i, int j) {
    T tmp = vec[i];
    vec[i] = vec[j];
    vec[j] = tmp;
}

template<class T>
void shuffleVector(vector<T> &vec) {
    for (int i = vec.size() - 1; i > 0; i--) {
        swap(vec, i, rand() % i);
    }
}

// This function converts decimal degrees to radians
double deg2rad(double deg) {
    return (deg * M_PI / 180);
}

//  This function converts radians to decimal degrees
double rad2deg(double rad) {
    return (rad * 180 / M_PI);
}

/**
 * Returns the earthDistance between two points on the Earth.
 * Direct translation from http://en.wikipedia.org/wiki/Haversine_formula
 * @param lat1d Latitude of the first point in degrees
 * @param lon1d Longitude of the first point in degrees
 * @param lat2d Latitude of the second point in degrees
 * @param lon2d Longitude of the second point in degrees
 * @return The earthDistance between the two points in meters
 */
double earthDistance(double lat1d, double lon1d, double lat2d, double lon2d) {
    double lat1r, lon1r, lat2r, lon2r, u, v;
    lat1r = deg2rad(lat1d);
    lon1r = deg2rad(lon1d);
    lat2r = deg2rad(lat2d);
    lon2r = deg2rad(lon2d);
    u = sin((lat2r - lat1r) / 2);
    v = sin((lon2r - lon1r) / 2);
    return 2.0 * EARTH_RADIUS * 1000 * asin(sqrt(u * u + cos(lat1r) * cos(lat2r) * v * v));
}

//double rad(double d)
//{
//    return d * PI /180.0;
//}
//
//double earthDistance(double lat1, double lng1, double lat2, double lng2)//lat1第一个点纬度,lng1第一个点经度,lat2第二个点纬度,lng2第二个点经度
//{
//    lat1 = 90 - lat1;
//    lat2 = 90 - lat2;
//    double a;
//    double b;
//    double radLat1 = rad(lat1);
//    double radLat2 = rad(lat2);
//    a = radLat1 - radLat2;
//    b = rad(lng1) - rad(lng2);
//    double s = 2 * asin(sqrt(pow(sin(a / 2), 2) + cos(radLat1) * cos(radLat2) * pow(sin(b / 2), 2)));
//    s = s * EARTH_RADIUS;
//    s = s * 1000;
//    return s;
//}

const gsl_rng *gslRng = gsl_rng_alloc(gsl_rng_taus);

double laplaceNoise(double epsilon) {
    // x ~ Laplace(0,1/epsilon)
    // https://www.johndcook.com/blog/2018/03/13/generating-laplace-random-variables/
//    double scale = 1 / epsilon;
//    double e1 = -scale * log(rand());
//    double e2 = -scale * log(rand());
//    return e1 - e2;
    return gsl_ran_laplace(gslRng, 1 / epsilon);
}

class Fingerprint {
public:
    int *values{};
    int size;

    friend ostream &operator<<(ostream &ostream, const Fingerprint &fp) {
        ostream << "(";
        for (int i = 0; i < fp.size; i++) {
            ostream << (to_string(fp.values[i]));
            ostream << ",";
        }
        ostream << ")";
        return ostream;
    }

    explicit Fingerprint(int prime_count) {
        this->size = prime_count;
        this->values = new int[this->size];
    }

    Fingerprint(const Fingerprint &fp) {
        this->size = fp.size;
        for (int i = 0; i < this->size; i++) {
            this->values[i] = fp.values[i];
        }
    }

    Fingerprint() : Fingerprint(PRIME_COUNT) {
    }

    ~Fingerprint() {
        delete[] values;
    }
};

struct fingerprint_equal {
    bool operator()(const Fingerprint &fp1, const Fingerprint &fp2) const {
        if (fp1.size != fp2.size) {
            return false;
        }
        for (int i = 0; i < fp1.size; i++) {
            if (fp1.values[i] != fp2.values[i]) {
                return false;
            }
        }
        return true;
    }
};

struct fingerprint_hash {
    size_t operator()(const Fingerprint &fp) const {
        unsigned long hash = 0;
        for (int i = 0; i < fp.size; i++) {
            hash = (fp.values[i]) xor hash;
        }
        return size_t(hash);
    }
};

vector<int> *generatePrimeList(int max) {
    bool *isPrime = new bool[max];
    auto primes = new vector<int>();
    memset(isPrime, true, max);
    for (int i = 2; i < max; i++) {
        for (int j = i; (long) i * j < max; j++) {
            isPrime[i * j] = false;
        }
    }
    for (int i = 1; i < max; i++) {
        if (isPrime[i]) {
            primes->push_back(i);
        }
    }
    return primes;
}

template<class T>
T &randomSelect(vector<T> *vec) {
    return vec->at(rand() % vec->size());
}

vector<int> *randomSelect(vector<int> *vec, int count) {
    auto result = new vector<int>;
    vector<int> copy(*vec);
    shuffleVector(copy);
    for (int i = 0; i < count; i++) {
        result->push_back(copy[i]);
    }
    return result;
}

double *millierConversion(double lat, double lon) {
    double L = 6381372 * M_PI * 2;
    double W = L;
    double H = L / 2;
    double mill = 2.3;
    double x = lon * M_PI / 180;
    double y = lat * M_PI / 180;
    y = 1.25 * log(tan(0.25 * M_PI + 0.4 * y));
    x = (W / 2) + (W / (2 * M_PI)) * x;
    y = (H / 2) - (H / (2 * mill)) * y;
    auto *result = new double[2];
    result[0] = int(x);
    result[1] = int(y);
    return result;
}

vector<string> *computeFingerprint(const char *str, int strLen, vector<int> &primeSuite) {
    int copyStrLen = strLen;
    const char *copyStr = str;
    auto fps = new vector<string>;
    for (int p : primeSuite) {
        strLen = copyStrLen;
        str = copyStr;
        int mod = 0;
        while (*str != 0 && strLen > 0) {
            mod = (mod * 256 + *str) % p;
            str++;
            strLen--;
        }
        fps->push_back(to_string(mod));
    }
    return fps;
}

string generateId(int len) {
    string id;
    for (int i = 0; i < len; i++) {
        int c = rand() % 62;
        if (c < 10) {
            id.append(to_string(c));
        } else if (c < 36) {
            stringstream ss;
            ss << static_cast<char>('a' + c - 10);
            id.append(ss.str());
        } else {
            stringstream ss;
            ss << static_cast<char>('A' + c - 36);
            id.append(ss.str());
        };
    }
    return id;
}

string generateId(int len, set<string> &notIn) {
    string id;
    do {
        id = generateId(len);
    } while (notIn.find(id) != notIn.end());
    return id;
}

set<string> *generateRandomIds(int count, int len, int seed) {
    srand(seed);
    auto idSet = new set<string>;
    string id;
    int c = 0;
    while (c < count) {
        id = generateId(len);
        if (idSet->find(id) == idSet->end()) {
            idSet->insert(id);
            c++;
        }
    }
    return idSet;
}

void blackhole(string str) {
    if (str == "abcdefghijklml") {
        cout << " ";
    }
}


// Simple RAII wrapper
class IntegrationWorkspace {
    gsl_integration_workspace *wsp;
public:
    IntegrationWorkspace(const size_t n = 1000) :
            wsp(gsl_integration_workspace_alloc(n)) {}

    ~IntegrationWorkspace() { gsl_integration_workspace_free(wsp); }

    operator gsl_integration_workspace *() { return wsp; }
};

// Build gsl_function from lambda
template<typename F>
class gsl_function_pp : public gsl_function {
    const F func;

    static double invoke(double x, void *params) {
        return static_cast<gsl_function_pp *>(params)->func(x);
    }

public:
    gsl_function_pp(const F &f) : func(f) {
        function = &gsl_function_pp::invoke; //inherited from gsl_function
        params = this;                     //inherited from gsl_function
    }

    operator gsl_function *() { return this; }
};

// Helper function for template construction
template<typename F>
gsl_function_pp<F> make_gsl_function(const F &func) {
    return gsl_function_pp<F>(func);
}

unordered_map<string, double> doubleIntegrationCache;
IntegrationWorkspace wsp1(1000);
IntegrationWorkspace wsp2(1000);

double integrateOverDoubleLaplace(double d, double ep1, double ep2) {
//    auto mxd = mxCreateDoubleScalar_800(d);
//    auto mxepi = mxCreateDoubleScalar_800(ep1);
//    auto mxepj = mxCreateDoubleScalar_800(ep2);
//    auto itg = new mxArray*[1];
//    itg[0] = mxCreateDoubleScalar_800(0);
//    mlfIntegrateDoubleLaplaceDistribution(1, itg, mxd, mxepi, mxepj);
    double epsabs = 1e-5;
    double epsrel = 1e-5;
    size_t limit = 50;
    double result, abserr, inner_result, inner_abserr;
    auto outer = make_gsl_function([&](double eta1) {
        auto inner = make_gsl_function([&](double eta2) {
            return ep1 * ep2 / 4 * exp(-ep1 * abs(eta1) - ep2 * abs(eta2));
        });
        gsl_integration_qagil(inner, eta1 - d, epsabs, epsrel, limit, wsp1,
                              &inner_result, &inner_abserr);
        return inner_result;
    });
    gsl_integration_qagi(outer, epsabs, epsrel, limit, wsp2, &result, &abserr);
    doubleIntegrationCache[to_string(ep1).append(to_string(d))] = result;
    return result;
//    return mxGetScalar_800(itg[0]);
}

int encrypt(string &clearText, unsigned char *key, unsigned char *ciphertext) {
    auto *iv = (unsigned char *) "0123456789012345";
    int len = 0, cipherLen = 0;
    EVP_CIPHER_CTX *ctx = EVP_CIPHER_CTX_new();
    EVP_EncryptInit_ex(ctx, EVP_aes_256_cbc(), nullptr, key, iv);
    EVP_EncryptUpdate(ctx, ciphertext, &len, reinterpret_cast<const unsigned char *>(clearText.c_str()),
                      clearText.size());
    cipherLen = len;
    EVP_EncryptFinal_ex(ctx, ciphertext + len, &len);
    cipherLen += len;
    EVP_CIPHER_CTX_free(ctx);
    return cipherLen;
}

string encrypt(string &clearText, string &key) {
    auto *cipherBuf = new unsigned char[128];
    int cipherLen = encrypt(clearText, (unsigned char *) key.c_str(), cipherBuf);
    return string((char *) cipherBuf, cipherLen);
}

int decrypt(string &cipher, unsigned char *key, unsigned char *clearText) {
    auto *iv = (unsigned char *) "0123456789012345";
    int len = 0, clearLen;
    EVP_CIPHER_CTX *ctx = EVP_CIPHER_CTX_new();
    EVP_DecryptInit_ex(ctx, EVP_aes_256_cbc(), nullptr, key, iv);
    EVP_DecryptUpdate(ctx, clearText, &len, reinterpret_cast<const unsigned char *>(cipher.c_str()),
                      cipher.size());
    clearLen = len;
    EVP_DecryptFinal_ex(ctx, clearText + len, &len);
    clearLen += len;
    EVP_CIPHER_CTX_free(ctx);
    return clearLen;
}

string decrypt(string &cipher, string &key) {
    auto *clearBuf = new unsigned char[128];
    int clearLen = decrypt(cipher, (unsigned char *) key.c_str(), clearBuf);
    return string((char *) clearBuf, clearLen);
}

string join(vector<string> *vec) {
    string str;
    for (const auto &item : *vec) {
        str.append(item);
    }
    return str;
}

template<class Key, class Value>
void memset(unordered_map<Key, Value> &map, set<Key> &keys, Value value) {
    for (const auto &k:keys) {
        map[k] = value;
    }
}

inline bool equals(const double &d1,const double &d2) {
    return abs(d1 - d2) < 1e-6;
}

#ifdef __MACH__

#include <mach/clock.h>
#include <mach/mach.h>

#endif

inline long long now() {
    timespec timestamp;
#ifdef __MACH__ // OS X does not have clock_gettime, use clock_get_time
    clock_serv_t cclock;
    mach_timespec_t mts;
    host_get_clock_service(mach_host_self(), CALENDAR_CLOCK, &cclock);
    clock_get_time(cclock, &mts);
    mach_port_deallocate(mach_task_self(), cclock);
    timestamp.tv_sec = mts.tv_sec;
    timestamp.tv_nsec = mts.tv_nsec;
#else
    clock_gettime(CLOCK_REALTIME, &timestamp);
#endif
    return (long long) timestamp.tv_sec * 1000000000LL + timestamp.tv_nsec;
}

#endif //LWPT_UTIL_H
